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Abstract: 



We determine the strong coupling a s at NNLO in perturbative QCD using the global 
dataset input to the NNPDF2.1 NNLO parton fit: data from neutral and charged current 
deep-inelastic scattering, Drell-Yan, vector boson production and inclusive jets. We find 
a s (Mz) = 0.1173±0.0007 stat , where the statistical uncertainty comes from the underlying 
data and uncertainties due to the analysis procedure are negligible. We show that the 
distribution of a s values preferred by different experiments in the global fit is statistically 
consistent, without need for rescaling uncertainties by a "tolerance" factor. We show 
that if deep-inelastic data only are used, the best-fit value of a s is somewhat lower, but 
consistent within one sigma with the global determination. We find that the shift between 
the NLO and NNLO values of a s is Aas° rt = 0.0018, and we estimate the uncertainty 
from higher-order corrections to be Aa^ NLO ~ 0.0009. 
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We have recently shown [T] that a surprisingly precise determination of a s using next- 
to-leading order perturbative QCD can be obtained from the dependence on the value of 
the strong coupling of the quality of a global unbiased fit [2] of parton distributions to 
deep-inelastic, Drell-Yan and inclusive jet data. This result was somewhat unexpected. 
Indeed, previous determinations of a s from global parton fits based on the idea of treating 
a s as a fit parameter [3H5] obtained larger uncertainties, and the precision of Ref. pQ is 
achieved even though the underlying PDF parametrization [2] adopted there is so general 
that the concept of a global best fit is no longer useful, and one must instead consider the 
fit quality (as measured by the x 2 ) as a Monte Carlo random variable which only stabilizes 
in the limit of very large Monte Carlo sample size. 

Given the great interest of an accurate determination of a s , in particular for Higgs 
searches at the LHC [6], it is important to repeat the analysis of Ref. [I] at NNLO, 
i.e. using the most accurate available QCD theory. This is the purpose of this study: 
besides determining a s and providing a detailed assessment of statistical and procedural 
uncertainties, and a general discussion of its reliability, we will provide an estimate of 
residual theoretical uncertainties. 

Our methodology is the same as in Ref. [T], to which we refer for a more detailed 
discussion. Here, it will be sufficient to recall that our determination is based on exploiting 
the fact that parton distributions determined with the NNPDF methodology [THE] are 
delivered for several values of a s , with all other aspects of the analysis being kept fixed. 
The quality of the fit of the data, as measured by its x 2 i can then be determined for each 
value of a s . Because NNPDF parton distributions are delivered as a set of Monte Carlo 
replicas, the x 2 1S a random variable, whose fluctuations only tend to zero in the limit of 
large sample size (as explicitly verified in Ref. [1]). The expected standard deviation of the 
X 2 for given sample size can be determined using the bootstrap method, as discussed in 
Ref. [1] . We thus arrive at a determination of the x 2 curve as a function of a s for a discrete 
set of values of a s and with an accuracy which depends on the size of the Monte Carlo 
replica sample. Both the set of values of a s and the sample size can then be enlarged, 
until satisfactory accuracy is reached. 

In Ref. [1] this was done by suitably enlarging both the number of values of a s and 
number of replicas of the available NNPDF2.1 NLO parton distribution sets [13]. Here, we 
start with the NNPDF2.1 NNLO parton distributions [2], which are available for a s (Mz) 
in the range from 0.114 to 0.124 in steps of 0.001, with N rcp = 1000 replicas for a s (Mz) = 
0.119, and -/V rcp = 100 replicas for all other values. We have extended this range down to 
0.110, and we have considerably increased the number of available replicas: we eventually 
end up with the number of replicas and values of a s summarized in Table HJ 

The x 2 values obtained for each value of a s are shown in Fig. [TJ the total number 



Values of a s (Mz) 




0.114,0.115,0.116,0.117,0.118,0.119 
0.112,0.113,0.120,0.121 
0.110,0.111,0.122,0.123,0.124 


1000 
500 
100 



Table 1: The number of PDF replicas N lcp in the PDF set for each value of a s (Mz). 
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Figure 1: The y 2 as a function of a s using NNLO NNPDF2.1 parton distributions. The uncer- 
tainty shown is due to the finite size of the PDF replica sample used to compute the x 2 f° r each 
value of a s . The curve is the result of a parabolic fit. 

of data points is N^ at = 3357 [2], and at the minimum x 2 /-^dat = 1-16. The error bar 
shown on each point is the standard deviation determined using the bootstrap method as 
described in Ref. pQ. A parabolic fit leads to 

nnlo ( Mz ) = Q.1173 ± 0.0007 stat ± 0.0001 proc , (1) 

where the statistical uncertainty corresponds to a shift by one from the minimum of the 
parabola, while the procedural uncertainty is propagated from the finite-size uncertainty 
on the x 2 values, and it is clearly completely negligible. The x 2 P er degree of freedom of 
the parabolic fit shown in Fig. [I] is Xpar/^dof = 1.1. This value remains unchanged if the 
outer two or four points are removed from the parabolic fit: this means that the behaviour 
of the curve shown is parabolic to within the accuracy of the individual points. 

Equation ([1]) is the main result of this paper. It is interesting to compare this to the 
NLO result of Ref. [I], namely 

c^ LO (M z ) = 0.1191 ± 0.0006 stat ± 0.0001 proc , (2) 

which was obtained using exactly the same methodology, from PDFs in turn determined 
using the same methodology and the same data 0; the value of the x 2 P er data point at 
the minimum at NLO is also the same, namely xV-^dat = 1-16. The result appears to 
be perturbatively quite stable. The small difference in statistical error between the NLO 

1 The number of datapoints for the NLO fit of Ref. [T] is iVdat = 3338 [13]: the slightly lower number of 
datapoints in the NLO case is due to the fact that the kinematic cuts on charm structure function data 
must be more stringent at NLO in order to exclude data for which NNLO corrections are very large, as 
discussed in more detail in Ref. [2]. 
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Experiment 




^.stat 


0-proc 


a tot 


P 


xV^dof 


NMCpd 


0.1146 


0.0105 


0.0019 


0.0107 


-0.25 


0.5 


NMCp 


0.1150 


0.0020 


0.0003 


0.0020 


-1.06 


1.6 


BCDMS 


0.1158 


0.0014 


0.0002 


0.0015 


-0.92 


1.2 


HERA-I 


0.1199 


0.0019 


0.0002 


0.0019 


1.31 


1.0 


ZEUS-H2 


0.1231 


0.0030 


0.0002 


0.0030 


1.89 


0.9 


NuTeV 


0.1177 


0.0038 


0.0004 


0.0039 


0.10 


1.2 


CDFZRAP 


0.1205 


0.0074 


0.0033 


0.0081 


0.40 


1.2 


CDFR2KT 


0.1225 


0.0021 


0.0003 


0.0021 


2.35 


1.8 


D0R2CON 


0.1111 


0.0028 


0.0004 


0.0029 


-2.10 


0.5 



Table 2: The best-fit values of a s obtained from a parabolic fit to the % 2 profiles of Fig. [2] Only 
experiments for which the parabola has a minimum in the fitted range are included. In each case 
we also show the statistical and procedural uncertainties (defined as in Eq. (H}), their sum in 
quadrature er tot , the pull P Eq. ([3]), and the x 2 of the parabolic fit. 

and NNLO fits indicates that the NLO fit quality deteriorates slightly faster when moving 
away from the minimum. However, this difference does not appear to be statistically 
significant: the gaussian uncertainty on the uncertainty is a a = a/y/2N, which in our 
case, with a = 0.0007 and N = A^of (the number of degrees of freedom of the parabolic 
fit) is (j (j = 0.0001. 

The statistical uncertainty on the determination of a s Eq. (fTJ) is very small when 
compared to the uncertainty on other existing determinations of a s : for example, the 
statistical uncertainty on the NNLO MSTW determination [3], which is also based on a 
PDF fit to a similar set of data, is twice as large. The size of the uncertainty on the 
MSTW08 result is significantly affected by the fact that a tolerance p3] criterion is used 
to obtain 68% confidence levels. This means that uncertainty ranges in parameter space 
are rescaled by an amount which is determined by requiring a reasonable behaviour of the 
distribution of best-fit results for the various parameters among the different experiments 
that enter the global fit, whose fluctuations may be significantly larger than expected on 
the basis of unrescaled uncertainties. It is thus important to establish whether such a 
rescaling is also necessary in our case. 

We do this (as in Ref. [I]) using a variant of the method suggested in Ref. [15]. Namely, 
we determine that the distribution of results obtained from individual datasets follows a 
Gaussian distribution. A possible deviation from unity of the width of this distribution 
would then indicate the need for tolerance rescaling. Of course care should be taken 
that a s values obtained from individual experiments in a global fit are correlated by the 
underlying global fit; these correlations could in principle lead to spurious minima in the 
X 2 profile which would be absent if a s was fitted from that experiment only, and thus 
would be of no significance. 

We have thus determined the profile of the contribution to the x 2 °f the global fit 
from each experiment entering the global NNPDF2.1 NNLO determination: all profiles 
are displayed in Fig. [21 along with their NLO counterparts of Ref. pQ. We have then 
determined a s and its uncertainty from a parabolic fit to each profile. For each experiment 
whose parabola has a minimum in the fitted range we have determined the pull, defined 
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Figure 2: Contribution to the total \ 2 from each experiment entering the NNPDF2.1 NNLO PDF 
determination, plotted as a function of a s , including the uncertainty on each point due to the finite 
size of the replica sample. The NNPDF2.1 NLO curves from Ref. [I] are also shown. 
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Figure 3: Histogram of the distribution of pulls included in Table [2j together with a gaussian fit 
to it. 



as 



Pi = 



af{M z ) - a s (M z ) 



(3) 



where a s (M z ) is the global best-fit value Eq. (P) and cr tot = J (a stat ) 2 + (cjP roc ) 2 the un- 



certainty on it determined from the statistical and procedural uncertainties also of Eq. (pQ), 
while a.g (M z ) and are the corresponding values for the i-th experiment. The best-fit 
value, statistical and procedural uncertainty for the i-th experiment are determined by 
applying to the x 2 profile for the contribution of the i-th experiment to the global fit 
exactly the same procedure used for the determination of the corresponding quantities in 
the global fit, and they are thus subject to the caveat that correlations with other experi- 
ments are neglected (for example, procedural uncertainties in each individual experiment 
will now have point-by-point correlation due to other experiments included in the global 



The best-fit values of a s , uncertainties and pulls are collected in Table along with 
the x 2 °f the parabolic fit to each individual profileU Clearly for some datasets (such as 
NMCpd or CDFZRAP) the minimum is actually quite shallow, as apparent from Figured 
and the uncertainty in Table[2]is correspondingly large: these data would not by themselves 
provide a determination of a s . However, note that here we are not determining the best-fit 
a s from each dataset, but rather studying the compatibility of different data sets within 
the global fit: from this point of view, data with a shallow minimum or no minimum are 
certainly compatible, and thus not problematic. 

The pull distribution has mean and standard deviation 



fit). 




a P = 1.46 ±0.24 



(4) 



2 Note that the definition of pull Eq. © is the same as in Eq. (3) of Ref. pQ, however all pulls were 
incorrectly given in Table 3 of that reference with the opposite sign convention. 



where the uncertainties on the mean and standard deviation are computed assuming that 
the population of nine pulls is extracted from a univariate Gaussian distribution. In Fig. [3] 
the histogram of the pull distribution is compared to a Gaussian normalized to the area 
under the histogram, and with mean and width given by Eq. We conclude that the 
pull distribution is consistent, with either no need or minimal need for tolerance rescaling. 
Even if one were to rescale the statistical uncertainty on a s by the pull, one would end up 
with a statistical uncertainty of 0.0010 instead of 0.0007. 

Having established that the estimate of the statistical uncertainty in our main result 
Eq. ([I]) is robust, we now turn to the issue of theoretical uncertainties. The main theoretical 
uncertainty comes from higher order QCD corrections. An indication of its size can be 
obtained by comparing results obtained at subsequent perturbative orders. Our NNLO 
and NLO results, Eqs. (HJ and ([2j) respectively, differ by 

A«P crt = K NL ° " «^ LO | = 0-0018. (5) 

This is likely to provide an upper bound to the uncertainty on the NNLO determination 
from higher order corrections. Hence, the perturbative stability of our result suggests that 
the higher order uncertainty on it is quite small. 

A somewhat less crude estimate can be obtained using a method recently proposed in 
Ref. [16], based on the idea that a Bayesian estimate of the behaviour of unknown higher 
orders of a perturbative series can be obtained from the behaviour of the known orders. 
The 68% confidence level (i.e. one sigma) uncertainty due to higher order corrections 
to the A:-th order result 

a k = cia l s + c m c4 +1 H h c k a k s (6) 

under the assumptions of the model in Ref. [16] is then given by 



a k s +l max{|ej|, . . . , |c fc |}^±±0.68 if 0.68 < ^ 

A *=< c (7) 

a k s +1 max{| Q |, . . . , | Cfc |} [(n c + l)0.32]~ 1/nc if 0.68 > 

with n c = k + 1 — I the number of known perturbative coefficients. 

In our case, the output value of a s in our analysis is not determined from the per- 
turbative expansion of a single quantity, but rather from the expansion for all processes 
which enter the global PDF fit. However, we can obtain an estimate for its perturbative 
behaviour by letting 

oo 

af(M z ) = J2c k [a° s (M z )} k (8) 

fc=0 

where af is the result of our determination at a given perturbative order, while a® is some 
underlying input "true" value of the strong coupling. We then know the values of at 
NLO and NNLO, but we do not know the value of at LO, nor the input value a®. We 
obtain our estimate by varying both in a wide range: we take a^ t,LO £ [0.100,0.180] and 
a° s £ [0.110, 0.125] and then we use Eq. ©, with af Eq. © evaluated at NLO or NNLO, 
and take in each case the largest output value for as the input parameters are varied. 
We then get 

A NLO = Ai = 0.0077 (9) 
a nnlo = A2 = 0.0009, (10) 
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Figure 4: Same as Fig. [TJ but using the NNPDF2.1 NNLO DIS only PDF set. 



so that in particular the NNLO uncertainty is about half the shift between the NLO and 
NNLO values. Clearly, results depend significantly on the range of variation of a^ t,LO : 
for example, if we restricted it to af' € [0.100,0.160] we would get Ai = 0.0053 and 
A2 = 0.0006. We have chosen a very broad range in order to obtain a conservative 
estimate. A more careful assessment could be obtained by performing a similar analysis 
(or factorization and renormalization scale variation) for each of the processes which enters 
the global fit, and then combining the ensuing uncertainties, while the analysis given here 
should be viewed as a semi-quantitative way of putting on a firmer footing the expectation 
that the uncertainty on a s at a given order due to higher-order corrections should be a 
fraction of the shift of the best-fit value of a s from the previous to the given order. 

We thus see that while at NLO theoretical uncertainties arising from higher perturba- 
tive orders are dominant over statistical uncertainties, at NNLO the two uncertainties are 
of similar size, and thus our determination of a s appears to remain quite accurate even 
when this uncertainty is included. Of the further sources of theoretical uncertainty, the 
largest is likely to be related to the treatment of heavy quarks, though such uncertainty 
is typically [T7] smaller than that related to higher order corrections. 

A sizable fraction of the data which are used for the NNPDF2.1 PDF determination 
comes from deep-inelastic scattering experiments. There is a widespread perception (see 
e -g- (18|J19|) that a s determinations based on deep-inelastic data lead to a value which is 
somewhat smaller than the global average. In Ref. [T], we have shown that in a DIS-only 
fit the BCDMS data (which have been long known [20J to favor a low a s ~ 0.114) have 
a minimum of the x 2 a t rather low a s < 0.113. However, this is no longer the case if 
the BCDMS data are included in a global fit along with Drell-Yan and jet data. We have 
traced this to the fact that at low a s DIS data and jet data pull the gluon in opposite 
directions, so the fit quality can be improved using DIS data only in a way which is 
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forbidden when jet data are present. Similar results have been obtained in Refs. [21(122]. 
where it was further argued that this effect is stronger if the gluon parametrization is too 
restrictive, and it was in fact shown explicitly that within the MSTW fitting framework 
low NNLO best-fit values (similar or lower to the values a s ~ 0.113 obtained at NNLO in 
Ref. [1]) are obtained if either the gluon parametrization is more restrictive, or jet data 
are excluded from the fit. In Ref. [1] we found that in the NNPDF framework (where an 
extremely flexible PDF parametrization is used) the effect of the BCDMS data on the DIS 
fit remains moderate: though the BCDMS data favor a low value of a s , the DIS-only best 
fit a s is still a s ~ 0.118, consistent with the global NLO best fit. However, the analysis 
of Ref. [1] is only performed at NLO, so one is immediately led to ask what happens at 
NNLO. 



NMCpd Rat NMCp SLAC 




0.109 0.111 0.113 0.115 0.117 0.119 0.121 0.123 0.125 0.109 0.111 0.113 0.115 0.117 0.119 0.121 0.123 0.125 0.109 0.111 0.113 0.115 0.117 0.119 0.121 0.123 0.125 
u,(M z ) o,(M z ) a,(M z ) 



Figure 5: Same as Fig. [2] but now comparing results obtained using the NNLO global and DIS- 
only fits. 

We have thus repeated our NNLO determination of a s , but now only including deep- 
inelastic scattering data in the computation of the x 2 - Starting with the NNPDF2.1 NNLO 
DIS-only fit of Ref. [2], we have generated N rep = 500 replicas for all values of a s (Mz) in 
the range from 0.110 to 0.124 in steps of 0.001. The x 2 values and finite-size uncertainties 
that we get are shown in Fig. HI to be compared to Fig. [T]of the global fit; the number of 
data points now is iVdat = 2783 [2], and at the minimum x 2 /-^dat = 114. A parabolic fit 
has a XparM'dof = 1.2, and leads to 

nnlo ( Mz ) = Q.1166 ± 0.0008 stat ± 0.0001 proc . (11) 
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NNLO a s ( M z ) from PDF Analyses 
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Figure 6: Comparison of the NNPDF2.1 NNLO a s (Mz) determinations to other recent determi- 
nations from NNLO PDF analysis. The PDG value of Ref. [T5] is also shown. 

Note that the outer two points at a s (Mz) = 0.110, 0.124 have not been used, because the 
quality of the parabolic fit would considerably deteriorate. The x 2 profiles from individual 
datasets in the global and DIS-only fit are compared in Fig. [5j 

Comparing our global and DIS-only determinations of a s Eqs. dH) and ([Til we see 
that, like at NLO [1], also at NNLO the DIS-only value, though somewhat smaller than 
the global one, agrees with it within the (small) statistical uncertainty at the one sigma 
level. The pattern of x 2 profiles for individual experiments is also similar at NLO and 
NNLO: the BCDMS data have a minimum for much lower a s in the DIS only fit than 
in the global fit. Note however that the overall DIS-only best-fit value of a s Eq. ([lip 
remains quite close to the global one, unlike in the MSTW08 case where the NNLO DIS- 
only value ends up being quite low (0.1139, to be compared to 0.1171 for the global fit) 
- perhaps because of our especially flexible PDF parametrization. In this respect, it is 
interesting to observe that MSTW find [21] that their DIS-only result depends strongly on 
whether the parametrization allows the gluon to become negative or not; in the NNPDF 
determination positivity is imposed at the level of physical observables, thus leading to a 
milder constraint on the PDF itself. 

Our results for the NNLO a s determinations, both from the global data set and from 
DIS data only, are displayed graphically in Fig. [6l where they are also compared to the 
other recent NNLO ABKM09 [4] and MSTW08 [3] determinations. Our results are in 
good agreement with MSTW, but with smaller uncertainties. The ABKM result is rather 
lower, but it should be kept in mind that ABKM only use DIS data and a smaller set of 
fixed-target Drell-Yan data, with a parton parametrization which is more restrictive than 
the MSTW one, so it is subject to the caveats mentioned above and discussed in Ref. [21]. 

In summary, we have presented the most accurate available NNLO a s determination 
from a global parton fit. The greater statistical accuracy of our result in comparison to 
other available determinations of a s along with PDFs [Hid] is due both to the use of the 
particularly flexible and reliable NNPDF approach to PDF parametrization and fitting, 
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and to the use of a wide and up-to-date dataset (in particular the recent combined HERA 
data [23] which are not used in Refs. [311]). The theoretical uncertainty on our determi- 
nation is likely to be dominated by higher-order perturbative QCD corrections, and thus 
there is no reason why it should be different from any other analogous determination (such 
as Refs. OH]). This uncertainty has not been determined accurately so far. However, the 
considerable perturbative stability we find, demonstrated by the moderate shift between 
our NLO and NNLO results, suggests that it might be smaller than hitherto expected: in 
particular, our estimate for it is less than half than that provided in Ref. [3]. It will be in- 
teresting to repeat this analysis once the very accurate data on standard candle processes 
from the LHC become available. 
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